% Compute statistics associated to the devaluation

Rn = SPn + RXn*wn; % sales in tradables
MFn= (1-sNn)*sMn * In; % final good imports in domestic labor
T =param.T;
disp('Check balanced trade condition after devaluation')
disp(((Mn+MFn)/(RXn+T)-1)*100)
disp('Growth in Wage and Spending')
disp([ (wn/w-1)*100 (Seq/SP-1)*100 ])

% Growth in total exports    
disp('Growth in total exports')
disp((RXn*wn-RX*w)/(RX*w)* 100)
%disp((TotalExports_norm1*weq-TotalExports_norm*w)/(TotalExports_norm*w)
%*100) % sanity check

Pn = (PriceChange/100+1)*P; % manufacturing price index in tradables
disp('% Growth in price of domestic input bundle')
disp(PriceChange)
disp('% Growth overall consumer price index')
disp((Ptotn/ (Ptot) - 1) * 100)

% change in correlation between sales and import shares conditional on
% importing
Dcorr1= corr(log(omega1(sD1<1)),1-sD1(sD1<1))-corr(log(omega(sD<1)),1-sD(sD<1));
   
disp('Change in aggregate import share and aggregate export share')
disp([(EqMoments.Aggsi-Aggsi)*100 (EqMoments.AggsX-AggsX)*100])   
disp('Change in correlation between sales and import shares')
disp(Dcorr1)
    
    
% Changes in distribution of import and export shares (Tables 8 and A2)
psi1 = prctile((1-sD1(sD1<1))*100,[10 25 50 75 90 95]); % percentiles of import shares conditional on importing
pxi1 = prctile(sX1(sX1>0)*100,    [10 25 50 75 90 95]);
dpsi =psi1-psi;dpxi =pxi1-pxi;
disp('Change in percentiles of import shares (importers)')
disp(dpsi([1 2 3 4 6]))
disp('Change in percentiles of export shares (exporters)')
disp(dpxi([1 2 3 4 6]))
    
si_avg1 = mean((1-sD1(sD1<1))*100); si_avg=  mean((1-sD(sD<1))*100); dsi_avg = si_avg1-si_avg;
xi_avg1 = mean((sX1(sX1>0))*100); xi_avg=  mean((sX(sX>0))*100); dxi_avg = xi_avg1-xi_avg;

% BHC Decomposition (Table 9)

[With, B, C, E, X, TOT,B2,Bex,Cex] = Decom(phi, sD, sD1,omega,omega1,sX1);
BHCdisplay = [With B Bex];
disp('BHC Decomposition: Within, Between, Between (Exporters)')
disp(BHCdisplay)
  
% Change in trade balance over sales (Table 10)
ExpoS = TotalExports_norm*w/R; % Exports normalized by Sales
ExpoS1 = TotalExports_norm1*wn/Rn;
ImpoS = (TotalImports_norm-MF)*w/R; 
ImpoS1 = (TotalImports_norm1-MFn)*wn/Rn;
disp('Change in trade balance over sales')
disp(((ExpoS1-ImpoS1)-(ExpoS-ImpoS))*100)

% Aggregate Elasticity of Substitution IMPORTS
disp ('Implied aggregate Armington elasticity')
growth_rs = ( EqMoments.Aggsi/(1-EqMoments.Aggsi)  / (Aggsi/(1-Aggsi)) - 1)*100;
EOS =1 - growth_rs / (-PriceChange);
disp(EOS)

% Check that labor in non-tradable sector is not in corner
if LN<param.L && LN>0
    disp('Interior allocation of labor in non-tradable')
else
    disp('Corner solution for allocation of labor in non-tradable')
end
if LNn<param.L && LNn>0
    disp('Interior allocation of labor in non-tradable (post devaluation)')
else
    disp('Corner solution for allocation of labor in non-tradable (post devaluation)')
end
    
% Extract results

micro = [omega,omega1,sD,sD1,sX,sX1];
Table5_first = [mean(log(fc)),mean((1-theta) * log(fx)), FM,FX,FXM, std(log(phi)),std(log(fc)),std(log(fx)), corr(log(phi),log(fc)),corr(log(phi),log(fx)), corr(log(fc),log(fx)), varrho, Timplied, param.theta]'; % Calibration Table
Table5_second = [Aggsi,AggsX,Firms_Importers+Firms_ImpExp,Firms_Exporters+Firms_ImpExp,Firms_ImpExp,STD_lnva,STD_si,STD_xi,corr_lnvasi,corr_lnvaxi,corr_sixi,sM,sN,(RXn*wn-RX*w)/(RX*w)*100]';
Table5_second = round(Table5_second,4);
Table5_first = round(Table5_first,3); 
Table6 = [psales;psi;pxi]; Table6 = round(Table6,2);% percentiles sales, si, sx
Table7 = [sicond;xicond;xicondsiXM]; Table7 = round(Table7,2); %Table 6 unconditional distributions
Table8 = [dpsi dsi_avg (EqMoments.Aggsi-Aggsi)*100 Dcorr1]; % changes in percentiles of import shares (importers), change in avg import share
Table8 = round(Table8,2);
Table9 = [round(BHCdisplay,2)];
Table10 = [(RXn*wn-RX*w)/(RX*w) (Mn*wn-M*w)/(M*w) ExpoS1-ExpoS-(ImpoS1-ImpoS)]*100; Table10 = round(Table10,2);
TableA2= [dpxi dxi_avg (EqMoments.AggsX-AggsX)*100]; TableA2 = round(TableA2,2); 
Text = [(wn/w-1)*100 (Ptotn/ (Ptot) - 1) * 100 PriceChange (Seq/SP-1)*100 EOS ]'; % growth in wage, in overall consumper price index (RER), in domestic input bundle price (P), and in domestic spending (S); also, aggregate armington elasticity
Text = round(Text,2);
Table_B5 = [(EqMoments.Aggsi-Aggsi)*100; Dcorr1; Table9(2); Text(5)];
    
cd ../../Output    
if glob==1
          csvwrite('../QuantitativeModel/temp/microdata_XM.csv',micro);
          csvwrite('Table5_first_XM.csv',Table5_first); 
          csvwrite('Table5_second_XM.csv',Table5_second); 
          csvwrite('Table6_XM.csv',Table6); 
          csvwrite('Table7_XM.csv',Table7); 
          csvwrite('Table8_XM.csv',Table8); 
          csvwrite('Table9_XM.csv',Table9); 
          csvwrite('Table10_XM.csv',Table10); 
          csvwrite('Stats_in_Text_Section4B_XM.csv',Text); 
          csvwrite('TableA2_XM.csv',TableA2); 
elseif ngf==1
          csvwrite('../QuantitativeModel/temp/microdata_NGF.csv',micro);
          csvwrite('TableA1_first_NGF.csv',Table5_first); 
          csvwrite('TableA1_second_NGF.csv',Table5_second); 
          csvwrite('Table8_NGF.csv',Table8); 
          csvwrite('Table9_NGF.csv',Table9); 
          csvwrite('Table10_NGF.csv',Table10); 
          csvwrite('Stats_in_Text_Section4B_NGF.csv',Text(5)); 
          csvwrite('TableA2_NGF.csv',TableA2); 
elseif monly==1
          csvwrite('../QuantitativeModel/temp/microdata_M.csv',micro);
          Table5_first(2) = 100000;
          csvwrite('TableA1_first_M.csv',Table5_first); 
          csvwrite('TableA1_second_M.csv',Table5_second); 
          csvwrite('Table8_M.csv',Table8); 
          csvwrite('Table9_M.csv',Table9); 
          csvwrite('Stats_in_Text_Section4B_M.csv',Text(5)); 
          csvwrite('TableA2_M.csv',TableA2);        
elseif nafta==1        
        csvwrite('TableA12_NAFTA.csv',Table8); 
        csvwrite('Footnote12_Appendix_NAFTA.csv',[Table9(1) Text(1)]);
elseif eta35glob ==1
          csvwrite('TableB5_Eta_PanelA_Col1.csv',Table_B5); 
elseif eta35NGF ==1
          csvwrite('TableB5_Eta_PanelB_Col1.csv',Table_B5); 
elseif eta045glob ==1
          csvwrite('TableB5_Eta_PanelA_Col3.csv',Table_B5);
elseif eta045NGF ==1
          csvwrite('TableB5_Eta_PanelB_Col3.csv',Table_B5);
elseif eps27glob ==1
         csvwrite('TableB5_Epsilon_PanelA_Col3.csv',Table_B5); 
elseif eps22glob ==1
         csvwrite('TableB5_Epsilon_PanelA_Col1.csv',Table_B5);
elseif eps27NGF ==1
        csvwrite('TableB5_Epsilon_PanelB_Col3.csv',Table_B5); 
elseif eps22NGF ==1
      csvwrite('TableB5_Epsilon_PanelB_Col1.csv',Table_B5); 
end
cd ../QuantitativeModel/mfiles